sa.dat.fn <- function(years,N0=1e6, n_intervals_year=12, M = 0.2, F = 0.3, q.s = 0.3, E.s = 0.05, ageplus = 40,
	sel.fn, catch.sel.pars = c(1,5), survey.sel.pars = c(5,2), samp.prob.c = 1, samp.prob.s = 1, N0.eq = 1e6, M0 =0.2, F0 = 0, pulse = FALSE){
	
	n_years <- length(years)
	catches <- surveys <- matrix(ncol = length(years)+1, nrow =12*length(years))
	colnames(catches) <- colnames(surveys) <- c(as.character(years), 'time.end')
	catches[,'time.end'] <- surveys[,'time.end'] <- years[1] + cumsum(rep(1/12,length(years)*12))
	
	P.eq <- cumprod(c(1,exp(-0.2*(1:(n_years-1)))))
	
	for(i in years){
		M.t <- rep(M,n_intervals_year*(years[n_years]+1-i))
		if(pulse) F.t <- rep(c(F,rep(0,11)),(years[n_years]+1-i))
		else F.t <- rep(F,n_intervals_year*(years[n_years]+1-i))
		E.s.t <- rep(c(rep(E.s,3),rep(0,9)),(years[n_years]+1-i))
		intrvls <- rep(1/n_intervals_year,n_intervals_year*(years[n_years]+1-i))
		x <- pop.sim.fn(N0=N0, A0=0, intervals = intrvls, F = F.t, M = M.t, q.s =q.s, E.s = E.s.t, 
			survey.sel.pars=survey.sel.pars, catch.sel.pars = catch.sel.pars, sel.fn = sel.fn)
		x <- cbind(x, time.end = i + cumsum(x[,'tau']))
		catches[,as.character(i)] <- x[match(catches[,'time.end'], x[,'time.end']), 'catch']
		surveys[,as.character(i)] <- x[match(surveys[,'time.end'], x[,'time.end']), 'survey']
	}
	
	N.eq <- N0.eq
	Z0 <- M + F0*sapply(0:(100-1) + 0.5, sel.fn, pars = catch.sel.pars)
	M.t <- rep(M,n_intervals_year*n_years)
	if(pulse) F.t <- rep(c(rep(F,1),rep(0,11)),n_years)
	else F.t <- rep(F,n_intervals_year*n_years)
	E.s <- rep(c(rep(E.s,3),rep(0,9)),n_years)
	intrvls <- rep(1/n_intervals_year,n_intervals_year*n_years)
	for(i in 1:100){
		N.eq[i+1] <- rbinom(1,size = N.eq[i], prob = exp(-Z0[i])) #stochastic eq. age comp?
		#N.eq[i+1] <- round(N.eq[i]*exp(-Z0[i])) #deterministic eq. age comp?
		x <- pop.sim.fn(N0=N.eq[i+1], A0=i, intervals = intrvls, F = F.t, M = M.t, q.s =q.s, E.s = E.s, 
			survey.sel.pars=survey.sel.pars, catch.sel.pars = catch.sel.pars, sel.fn = sel.fn)
		x <- cbind(x, time.end = years[1] + cumsum(x[,'tau']), age.mid = i + cumsum(intrvls)-intrvls/2, cohort = years[1] - i)
		if(i == 1) {
			catches.eq.pop <- x[which(x[,'catch']>0),c('time.end','age.mid','catch','cohort')]
			surveys.eq.pop <- x[which(x[,'survey']>0),c('time.end','age.mid','survey','cohort')]
		}
		else{
			catches.eq.pop <- rbind(catches.eq.pop, x[which(x[,'catch']>0),c('time.end','age.mid','catch','cohort')])
			surveys.eq.pop <- rbind(surveys.eq.pop, x[which(x[,'survey']>0),c('time.end','age.mid','survey','cohort')])
		}
	}
	
	temp <- cbind(time.end = catches[,'time.end'], age.mid = catches[,'time.end'] - years[1] - 0.5/n_intervals_year, n.catch = catches[,1], cohort = years[1])
	temp <- temp[which(catches[,1]>0),]
	catches.rf <- temp
	for(i in 2:n_years){
		temp <- cbind(time.end = catches[,'time.end'], age.mid = catches[,'time.end'] - years[i] - 0.5/n_intervals_year, n.catch = catches[,i], cohort = years[i])
		temp <- temp[which(catches[,i]>0),]
		catches.rf <- rbind(catches.rf, temp)
	}
	print(dim(catches.rf))
	print(head(catches.rf))
	catch.raw <- rbind(catches.rf, catches.eq.pop)
	if(samp.prob.c<1) catch.raw <- cbind(catch.raw, n.aged =rbinom(dim(catch.raw)[1],size = catch.raw[,'n.catch'], prob = samp.prob.c))
	else catch.raw <- cbind(catch.raw, n.aged = catch.raw[,'n.catch'])
	temp <- catch.raw
	y <- round(temp[,'time.end'],4)
	temp[which(y %in% c(years,years[n_years]+1)),'time.end'] <- temp[which(y %in% c(years,years[n_years]+1)),'time.end'] - 0.0001
	caa <- tapply(temp[,'n.aged'],list(floor(temp[,'time.end']), floor(temp[,'age.mid'])), sum)
	caa[is.na(caa)] <- 0
	ind <- which(as.integer(colnames(caa))>= ageplus)
	caa <- cbind(caa[,-ind], apply(cbind(caa[,ind]),1, sum))
	colnames(caa)[dim(caa)[2]] <- ageplus
	
	catch.y <- tapply(temp[,'n.catch'],list(floor(temp[,'time.end'])), sum)
	
	temp <- cbind(time.end = surveys[,'time.end'], age.mid = surveys[,'time.end'] - years[1] - 0.5/n_intervals_year, n.catch = surveys[,1], cohort = years[1])
	temp <- temp[which(surveys[,1]>0),]
	print(temp)
	surveys.rf <- temp
	for(i in 2:n_years){
		temp <- cbind(time.end = surveys[,'time.end'], age.mid = surveys[,'time.end'] - years[i] - 0.5/n_intervals_year, n.catch = surveys[,i], cohort = years[i])
		temp <- temp[which(surveys[,i]>0),]
		surveys.rf <- rbind(surveys.rf, temp)
	}
	print(dim(surveys.rf))
	survey.raw <- rbind(surveys.rf, surveys.eq.pop)
	if(samp.prob.s<1) survey.raw <- cbind(survey.raw, n.aged = rbinom(dim(survey.raw)[1], size = survey.raw[,'n.catch'], prob = samp.prob.c))
	else survey.raw <- cbind(survey.raw, n.aged = survey.raw[,'n.catch'])
	temp <- survey.raw
	y <- round(temp[,'time.end'],4)
	temp[which(y %in% c(years,years[n_years]+1)),'time.end'] <- temp[which(y %in% c(years,years[n_years]+1)),'time.end'] - 0.0001
	saa <- tapply(temp[,'n.aged'],list(floor(temp[,'time.end']), floor(temp[,'age.mid'])), sum)
	saa[is.na(saa)] <- 0
	ind <- which(as.integer(colnames(saa))>= ageplus)
	saa <- cbind(saa[,-ind], apply(cbind(saa[,ind]),1, sum))
	colnames(saa)[dim(saa)[2]] <- ageplus
	
	survey.y <- tapply(temp[,'n.catch'],list(floor(temp[,'time.end'])), sum)

	
	return(list(catches = catches, surveys = surveys, N.eq = N.eq, catches.rf=catches.rf, surveys.rf=surveys.rf, catches.eq.pop = catches.eq.pop, 
		surveys.eq.pop = surveys.eq.pop, survey.raw = survey.raw, catch.raw = catch.raw, saa = saa, caa = caa, survey.y = survey.y, catch.y = catch.y))
}
